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Abstract. When submillimetric particles saturated with a fluid form a 
compact cluster on the top on a confined clear fluid of the same na- 
ture, they fall by detaching from the lower boundary of the cluster, and 
this separation front between particles and fluid is unstable. Particles 
detach and fall in the clear fluid region, giving rise to growing fingers 
of falling particles. We study this problem using both experiments and 
hybrid granular/fluid mechanics models. In the case of particles from 50 
to 500 microns in diameter falling in air, we study the horizontal den- 
sity fluctuations at early times: the amplitude of the density difference 
between two points at a certain horizontal distance grows as a power 
law of time. This happens up to a saturation corresponding to a power 
law of the distance. The way in which the correlation length builds up 
to this saturation also follows a power law of time. We show that these 
decompaction fronts in sedimentation problems follow a Family- Vicsek 
scaling, characterize the dynamic and Hurst exponent of the lateral 
density fluctuations, respectively z ^ 1 and ( ^ 0.75, and show how 
the pref actors depend on the grain diameter. We also show from simi- 
lar simulations with a more viscous and incompressible fluid, that this 
feature is independent of the fluid compressibility or viscosity, ranging 
from air to water /glycerol mixtures. 



1 Introduction 



In many natural and industrial situations, mixtures of granular materials and fluids 
are deformed. Physicists have intensively studied granular mechanics using discrete 
methods over the last 30 years [1,2^ since the rise of powerful computers, with various 
types of avalanche dynamics [3], the dynamics of sheared layers |4|5j and the stress 
distribution in sandpiles [6^. For sufficiently small grains, when the fluid seeps through 
the granular assembly, the resulting drag is sufficient to rearrange the grains, and it 
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is necessary to take into account the momentum exchange between the granular and 
the fluid phase to understand the dynamics. 

Such hybrid granular/fluid flows are important /manifest in many natural flows at 
the surface of the earth: e.g. in liquefaction [7 , mudfiows |8|9|lQj . or sand volcanoes. 
They also play a role in underwater avalanches and turbidites, or in sedimentation 
processes when large amounts of sediments are released and deposited. Eventually, 
their understanding is/are crucial for the mechanical stability of/in borehole exploita- 
tion/exploration, during pumping in to or out of man-made wells, in hydrofracture 
[11], or for the understanding of reservoir stability and fault lubrication [12 . 

Elementary situations of sedimenting grains/unconsolidated grains (in a fluid) 
and their stability have been recently studied in generic situations where a dense 
granular packing saturated with a fluid is subjected to a flow of the same fluid. 
For example, the propagation of rising bubbles in air- grain flows in tubes, or the 
intermittency of flows in ticking hourglasses has been studied [1 3|14|15j . Experiments 
with dense granular packings displaced by pressurized air have been followed optically 
in Hele-Shaw cells [1 6|17|18|19|2Qj . and studied numerically using models where the 
momentum exchange is based on a Darcy description, with a flow described at a 
coarse grained scale of a few grain diameters [ 21|22j . 

For example, when a granular packing in a circular cell where gravity acts per- 
pendicularly to the cell walls starts to flow due to a sufficient air-overpressure in 
the centre, an instability corresponding to a granular Saffman- Taylor analog was ob- 
served |16|23 j. Such analog behavior was also observed in linear cells, either with a 
free boundary at the outlet leading to decompaction |17|18j , or with a confined flow 
for the particles leading to compaction |24|25j. Changing the properties of the inter- 
stitial fluid with orders of magnitudes in viscosity and compressibility (between oil 
and air) was shown to affect the dynamics of this process, but to lead to the formation 
of similar patterns \T8* . 

When gravity is acting along the flow, a granular analog of the Ray leigh- Taylor 
instability has been studied: situations where initially a pack of dense grains are re- 
leased into a zone of lighter fluid. The same fluid is also present in the interstitial 
space between the packed grains. This situation would naturally arise during sedi- 
mentation flows, e.g. when a sudden amount of compact material is released in open 
water. This instability has been studied with mixtures of grains (typically 100 pm in 
diameter) and air |19|26|27j , and with less compressible and more viscous fluids such 
as water glycerol or oil [2 0l28l29j . 

The release of such dense packs of grains in a clear fluid contained in a closed cell, 
leads to the formation of a detachment front at the bottom of the pack. The particles 
fall faster from the bottom part than from the upper part of the pack, and fingers form 
in the previously clear fluid region. This front is unstable, grows fast in amplitude, 
and leads to coarsening bubbles empty of grains rising through the grain pack. The 
dynamics of this front has been studied using Fourier analysis by investigating the 
growth of each mode [19^ 26 30 31 . In this system, the change of the particle size, for 
sizes ranging from roughly 50 to 500 /im in diameter, was shown to leave the system 
dynamics invariant, using the grain diameter as the characteristic length scale. 

We will show in this manuscript that the horizontal density profiles, as function 
of time and space, for different particle sizes, can be entirely collapsed and follow a 
Family- Vicsek scaling. Such scaling, introduced by Family and Vicsek [32,, was found 
to be obeyed by many surface growth processes. It relates the anomalous diffusion of 
a quantity, the growth of correlations, and the self affine geometry of the developed 
interface - or of the quantity characterized by this scaling [33 . For example, such 
a collapse was observed for the thermal fluctuations of magnetic hole assemblies, 
which are dipolar chains in thin fluid cells j 34|35|36] . It is also obeyed by natural 
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Fig. 1. (a) Experimental setup and principle, (b) typical granular Ray leigh- Taylor flow, 
where fingers of falling grains (white) form in the previously air-filled region (black) during 
an experiment, and (c) the experimental cells of different sizes with grains of proportional 
diameters. 

patterns slowly growing in rocks, as pressure solution seams evolving into stylolites 
|37l38l39l40| . 

This implies that a correlation length grows in this system, as a power law of time, 
^ ~ t^/^, with a dynamic exponent around z ^ 1. For scales above ^, the root mean 
square of the lateral density fluctuations grows in a super diffusive way, as a ~ 
with a growth exponent around P ~ 0.75. Below the correlation length, the density 
profiles are shown to be self-affine, with a Hurst exponent around ( ~ 0.75. 

In the next section, we will summarize the methodologies used to study the gran- 
ular Rayleigh- Taylor instability, from an experimental and from a numerical perspec- 
tive. Next, we will describe the types of patterns obtained at different particle sizes. 
Eventually, we will perform the Fourier analysis of the density fluctuations, and show 
that they can be collapsed at different times, for different wave numbers and for all 
investigated particle sizes on a master curve, displaying the fact that they follow a 
Family- Vicsek scaling. 



2 Rayleigh-Taylor instability: definition and methods 
2.1 Experimental setup 

Historically, the Rayleigh-Taylor instability arises from a flat interface when a dense 
fluid is initially situated above a lighter one: in a closed container, the two fluids 
interchange positions while forming a fingering pattern. This instability, and the se- 
lection of a preferred wavelength due to the interplay between gravitational, capillary 
and viscous forces was first studied by Lord Rayleigh and L G. Taylor |41|42j . The 
granular analog of this instability has been studied recently with air or more viscous 
fluids: it arises when initially a pack of dense grains is released on top of clear fluid 
in a closed vessel. There is no surface tension in this problem since the fluid is also 
initially present in the porous space between the grains in the upper pack. 

The experiments are performed in a closed impermeable Hele-Shaw cell, which 
is initially prepared by letting grains accumulate at the bottom. The cell is then 
suddenly rotated to be brought upside down, and a damping stops the cell in a 
vertical direction, attenuating vibrations as much as possible. The setup is illustrated 
in Fig. [l] Several system sizes have been used for air-grain experiments, ranging from 
80 /i m diameter beads to 570 /im diameter ones for the air/grain system. The system 
sizes are proportional to the bead diameters, from 31 to 200 mm. The grains used 
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Fig. 2. Experimental snapshots of the granular Ray leigh- Taylor instability in air where the 
four cells shown in Fig. [T]^c) are used. 



in the experiments where polystyrene spheres (Ugelstad spheres from Microbeads), 
with a density of 1.05 g/cm^. 

Early stages of the experimental pictures are shown on Figs, [l] (right) and Fig. |2] 



2.2 Simulation technique 



The simulation techniques are described in details in [19* for the flow of grains and 
compressible air, or [20 for the flow of grains and incompressible viscous fluid. 

In summary, the principle is as follows: the grains interact via repelling contact 
forces and via a drag force due to the flow of the interstitial fluid. This can be 
expressed as Newton's second law for each grain, as 

ma = Fj + {pg - pf)g - VP/p (1) 

where m and a are the grain mass and acceleration, Fi corresponds to the solid 
interaction forces with neighboring grains, the third term corresponds to buoyancy 
forces, and the last one to seepage forces: g is the gravity, pf and pg are the fluid 
and solid bulk mass densities, P = Pt — pjgz is the fluid pressure deviating from the 
hydrostatic profile - where Pt is the total fluid pressure, p/ is the fluid density, g 
the gravity and z the depth. This pressure is evaluated over a grid discretized at the 
representative elementary volume scale (the Darcy scale, a few grains large), and p 
is the number density of the grains: p = Vg/{1 — ^) with Vg the grain volume, and <P 
the porosity. 

This equation is being solved at each time step, using a Verlet algorithm in molec- 
ular dynamic codes. The forces Fj are determined from elastic contact model [16 20 , 
or from a contact dynamics scheme [19]. 

The dynamic equation solved for the fluid results from the following principles: 
mass conservation for the grains, mass conservation for the fluid, equation of state 
for the fluid, and Darcy to get the relative velocity between the fluid and the solid. 
Calling u the grain velocity, v the fluid one, and pf the fluid density, this can be 
expressed at the Darcy scale as: 

S^(l-^) + V-((l-^)ix) = 
dt{pf^)^V - i^pfv) = 

K.\/P 

V — u = — 

p<p 



(2) 
(3) 

(4) 
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where /i is the fluid viscosity, is the local permeability of the packing, evaluated 
from the Carman-Kozeny equation: 

180(1-^)2 

and the fluid mass conservation is expressed as one of the two alternatives: in general, 
linearizing this equation of state, one gets 

Pf{P)=Po{l + /3{P-Po)) 

where /3 is the fluid compressibility, Pq an initial (reference) pressure, and po the fluid 
density at that pressure. 

This set of equations leads to the following reduced form: 

^ {dtP + u • VP) = V • (^P^\/P^ -PV 'U (5) 

where P = p//(po/5) = P - Pq + 

Considering air as a perfect gas, this reduces to = poP/Po, = 1/Po and 
P = P, so that 

^ {dtP + u • VP) = V • (^P^\/P^ - PV -u (6) 

which is solved alternately with the molecular dynamic time steps, using a Cranck- 
Nicholson algorithm [43 . 

For a viscous fluid considered as incompressible, we have /3 = 0, p/ = po- and the 
equation [5] in the limit P ^ 00, reduces to 

V- (^^VP^ = V -u 

which is solved at each time step using a multigrid algorithm [43]. 



2.3 Dynamics of the granular Rayleigh-Taylor instability 

In both experiments and simulations, for grains falling in air or in water/glycerol, 
the dynamics can be qualitatively described as follows: From the initial configuration 
where grains are packed above the clear fluid region, the dense grains start to flow 
downward, i.e. sediment, while the fluid passes through the grains. This situation 
leads to the formation of fingers of particles detaching from the lower boundary of the 
upper compact pack. These fingers of particles regroup, and bubbles of low particle 
density form in the pack, rise and coalesce. This dynamics is illustrated in Fig. [2] 
(experimentally) and Fig. |3] (numerically), for the air/grain case, and in Fig. [4] for the 
incompressible fluid/grain case. 

We notice that at low viscosities, simulations with incompressible fluids are remi- 
niscent of the dynamics obtained with compressible air. As the viscosity is increased 
(from top to bottom in Fig. |4|, the fingers get mixed more rapidly as they detach 
from the lower front of the compact pack. The grains have here a size of 140 microns 
in diameter, and a density of 2.5 g/cm^, the fluid has the density of water. 

One way to analyze the instability of the detachment front, at the bottom of the 
compact pack, is to compute horizontal density profiles. In this way, it is possible 
to measure the lateral growth of density fluctuations, while the fingers form and 
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Fig. 3. Simulation at four different grain/system sizes (horizontal axis) at two different times 
(vertical axis) showing the development of the instability with the initial front destabilization 
in the top row. 



regroup. This can be done at each horizontal level at vertical coordinate y and time 
t, extracting density fronts as function of the horizontal variable pix^y^t). This is 
illustrated in Fig. |5j 

The most variable profiles as function of x are just below the detachment front: 
the particle density is roughly constant (at a scale above a few grains) in the compact 
packing above the detachment front, and constant and low in the empty region below. 
Hence, to increase the statistical quality of the quantity we analyze, we perform the 
average over all horizontal profiles at a given time, and obtain in this way horizontal 
density profiles p{x^t) = {p{x^y^t))y. For the experiments, the density is not known 
directly, but the same procedure is performed on the basis of light intensity (pixel 
gray value) recorded with a fast camera. 

These horizontal density functions are then spatially Fourier transformed, and 
their Fourier power-spectra are represented at various times and for various grain 
sizes using a bi-logarithmic scale, as is shown in Fig. [6] for the simulations, and in 
Fig. [7| for the experiments. The inset plots represent the raw power spectra for all 
sizes, and the full figures indicate the power spectra using the grain diameters as a 
unit size. The comparison between both shows the size collapse. 

The experiments and simulations with compressible air have been analyzed as 
function of the particle size. Snapshots of the early times of these simulations with 
variable particle sizes are shown in Fig. [2] A spatial rescaling using the grain diameter 
rendered the pattern formation mechanism similar for the range of diameters investi- 
gated. This was studied in details in [27 . The rescaling is illustrated by the collapse 
of the structures analyzed in Fig. [g] (for simulations) and Fig. [t] (for experiments), 
obtained after rescaling the spatial units with the grain diameter, d. The insets in 
these figures show the density profiles prior to rescaling. 
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Fig. 4. Numerical simulations of grains sedimenting in an incompressible fluid of variable 
viscosities /i/ (increasing from top to bottom). Time runs from left to right. 
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wave number {^/cwi 



Fig. 5. Extraction of a density function under the dense plug (middle dashed line in the 
right figure), and computation of its Fourier power spectrum (left). 



10^ 



10 



10" 



10 











a) 

. t = 0.14 s 






























10° 








10"^ 












10 ' 


lo" 


lo' 



10^ 



10 



10" 



10 



10"^ 



dk 



10" 




70 |im 

140 |im 

210 |im 
280 |im 
350 |im 
420 |im 
- - - 490 |im 



Fig. 6. Power spectrum collapse, simulations: spatial scaling. 
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Fig. 7. Power spectrum collapse, experiments: spatial scaling. 

3 Family- Vicsek scaling of the lower decompaction front 

In the plots of the Fourier transformed density profiles in Figs. |6] and [7j one notes 
that, both for the experiments and for the simulations, the low wavenumber part of 
the spectra is flat. This corresponds to spatially uncorrelated density fluctuations for 
large scales (with a white noise character), whereas the power spectrum decreases as 
a power-law for small scales (high wavenumber). A behavior 
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is observed on this part, both for experiments and for simulations. The cutoff kc 
corresponds to the maximum size of the density fluctuations, ^ = 27r /kc^ called the 
correlation length. The power-law behavior for scales smaller than the correlation 
length corresponds to the self-affine character of these density fluctuations. The power 
spectrum distribution of a self-affine quantity may be expressed as 

S{k,t)^k-^'^^^\ 

which in our case yields a Hurst exponent ( ^ 0.75. 

Deviations are observed for the experiments in Fig.[7|on the smallest scales, which 
can be attributed to a different dynamics at scales of the order of the plate separation. 

Interestingly, when the case of incompressible fluids is compared, a similar be- 
havior is observed as for the compressible case, see Fig. [Sj This behavior extends to 
the whole range of viscosities probed, from 0.018 to 0.9 mPa-s, almost three orders 
of magnitude. This is the case despite the obvious difference between the developed 
stages of the instability, as is seen by comparing Figs. |3] and [4j The early stages, how- 
ever, seem to follow the same type of behavior with respect to the lateral spreading 
of the density fluctuations. 

The Family- Vicsek scaling behavior is summarized in the plots in Fig. [9] for the 
compressible case, and in Fig. [lO] for the incompressible case. For the compressible 
case, the transformed density functions S{k^t^d) for all the probed grain sizes d and 
for times t < 0.23 s all fall on a single master curve shown in Fig. |9|a). For the 
incompressible case there is, however, a deviation from the Family- Vicsek scaling for 
the early times and largest viscosities probed, as seen in the incomplete collapse for 
the first times in Fig. 10 'c). 



The observed Family- Vicsek scaling merely reflects the following about the power 
spectrum S{k^ t): 

S{k, t) = Ck-^-^^g{dk{t/r)^^') (7) 
g{x) - Cix^^^^ for x <Cl, i.e. l//c> ^ with ^ - t^^^ (8) 
g{x) = C2 for X > 1, i.e. 1/fe <C ^ , (9) 

where r is a time-constant and and C2 are constants. So that S{k^ t) is inde- 

pendent of k for scales above the correlation length ^, and independent of t when the 
mode k has saturated, with a wavelength 1/k exceeding the correlation length. This 
leads to a generalized diffusion behavior 

S{k,t)^t^^^^^^/' (10) 

when 1/k > ^, and 

S{k,t) ^ k-^^^^^^ (11) 

when 1/k < ^. 

This corresponds, in real space, to a behavior for the root mean square of the solid 
density 

following a Family- Vicsek scaling that can be obtained form Parseval's theorem, as 



POO 

a^{L,t) = (1/L) / S{k,t)dk 

Jk=2iT/L 



(12) 



= (1/L) / Ck-^-^^g{dk{t/ry/')dk (13) 

Jk=27v/L 



00 



(t/T)(2C)A(i/L) / k'g{k')dk', (14) 

J27r(t/r)i/^/L 
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Fig. 8. Power spectra of the horizontal density profile for different times: results of simula- 
tions in an incompressible fluid. The viscosities of the fluids are respectively (a): 0.018, (b): 
0.073 and (c): 0.9 mPa-s. The different curves correspond to different times after the start 
of the instability. 




(dk)*t [s] (dk)*t [s] 



Fig. 9. Family- Vicsek scaling obtained from numerical data, (a) shows the scaling for all 
grain sizes and for t < 0.23 s, while (b) is the same plot as (a) but only showing the 210 pm 
data and where the colors indicate the time from blue (early times) to red (late times). 
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amt) = it/rfF{L/t'/^) 

with 

l3 = C/z. 

Hence, the density profile under the detachment front saturates as a self-afRne func- 
tion [33 with a Hurst exponent ( 0.75 for scales below ^, and displays a superdiffu- 
sive behavior with an exponent /3 :^ 0.75. The root mean square of the solid density, 
— (P^ ~ (p)^)^^^ displays a behavior corresponding to 

Conversely, for scales / above ^, one has a self-affine behavior, 

{{p{x + l)-p{x)f^ei. 

The correlation length ^ is increasing as a power law with time, as ^ ^ t^/^. The 
density in the detachment front exhibits these behaviors with respectively dynamic. 
Hurst and growth exponents 

z = 1 (15) 
C = 0.75 (16) 
(3 = C/z = 0.75 (17) 

One can also equivalently express the Family- Vicsek scaling with a scaling function 
such as: 

S{k, t) = C{dk)-^-^g{dk{t/r)) (18) 
g{x) - Cix^-^ forx <C 1, i.e. 1/k > ^ with ^ - t (19) 
g{x) = C2 forx > 1, i.e. 1/A: <C ^. (20) 

Indeed, S{k^t) * (d* /c)^*^ = Cg{d^ k^t/r) is represented for air/grain simulations in 
Fig. [9) and this scaling function master curve displays these two expected behaviors 
for the granular-gas flow. 



This is also the case for the granular-incompressible fluid flows, as shown in Fig. 10 
for the three viscosities. As observed in the previous representation, the collapse is 
rather followed, apart from the earliest times at larger viscosities. 



4 Conclusion 

We have studied situations where well packed cohesionless grain assemblies, heavier 
than the fluid between them, are released in a clear fluid. This gives rise to the so- 
called granular Ray leigh- Taylor instability, where bubbles empty of grains rise, while 
granular fingers form in the clear fluid region. 

The principles of hybrid simulations starting from basic physical equations was 
established, for both cases where the carrier fluid is a perfect gas or an incompressible 
viscous fluid. The experimental setup and results are also shown to be consistent with 
the simulations in the air /grain case. 

Grains are mixed during this process, and one can follow the evolution of spatial 
correlations in lateral density fluctuations, with vertically averaged horizontal grain 
density profiles p{x). These density fluctuations are shown to follow an anomalous 
diffusion behavior with a growth exponent P ~ 0.75. Hence, for two points at a large 
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Fig. 10. Master function intervening in the Family- Vicsek scaling for the Ray leigh- Taylor 
instability in incompressible fluids: Scaled power spectra S{k, t)*((i*/c)^"^, which are functions 
gofd^k^t/r, with a behavior g(x) ^ Cix^'^ at small arguments, and g{x) ^ C2 at large 
ones. The Hurst exponent used for the scaling is 0.75, and the dynamic exponent is 1. The 
viscosities of the fluids are respectively (a): 0.018, (b): 0.073 and (c): 0.9 mPa-s. 



horizontal distance x from each other, a difference in density will grow as 5p{x^ t) ~ 
. The density fluctuations will follow this behavior up to a saturation when the 
correlation length ^ reaches x. ^ is growing as a power-law of time, with a dynamic 
exponent around 1. Eventually, for scales larger than ^, the density profile displays 
a self-affine behavior, with a Hurst exponent around 0.75. This is seen in the spatial 
Fourier domain, both for experiments and simulations for air-grain systems, and for 
simulations for mixtures of grains and incompressible viscous fluid. 

Eventually, we have shown that the power-spectra of the density at all times can 
be collapsed according to a Family- Vicsek scaling. 
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